******************************** Set Dirs ***************************************


 
use .\data\clean\revision1\5minpemsspeedflowforRD,clear
set matsize 5000
 ******************************** Prep Data ***************************************

drop if inlist(dow,6,0)                  // No weekends
drop if holiday==1 				// No holidays
sort date hour
merge m:1 date hour using ".\data\clean\revision1\\gprice.dta", keep(matched) nogen
sort station date hour 
merge m:1 station date hour using ".\data\clean\revision1\\lagged_210_hourlyspeed.dta", keep(matched) nogen
sort station date hour 
merge m:1 station date hour using ".\data\clean\revision1\\weather.dta", keep(matched) nogen


gen ln_flow = log(flow)
gen ln_speed= log(speed)
gen ln_210  = log(speed210)
gen ln_gp   = log(gprice)
gen vis2    = vis^2
gen prec2   = prec^2

drop if date<mdy(12,01,2012) /*Opening of second HOV lane */
drop if date>mdy(12,30,2013)
g period = cond(hour>4 & hour<9,"AM",cond(hour>8 & hour<16, "OFF1", cond(hour>15 & hour<20, "PM", "OFF2")))
gen EL      = date>=mdy(02,23,2013)

*Cluster Variable Week Count
g year=year(date)
g month=month(date)
gen weekid = week(date)
replace weekid=weekid+100 if year==2013
replace weekid=weekid+200 if year==2014
replace weekid=weekid-100 if year==2011

egen t=group(date)


qui sum t if date==td(22feb2013)
g poly=t-`r(max)'-1

***************************************** Local Linear Regressions **************************************************
gen w55=max(1-(abs((date-mdy(02,23,2013))/55)),0) // Triangular Kernel
gen w60=max(1-(abs((date-mdy(02,23,2013))/60)),0) // Triangular Kernel
gen w65=max(1-(abs((date-mdy(02,23,2013))/65)),0) // Triangular Kernel
*gen ew=max((3/4)*(1-(((date-mdy(02,23,2013))/60)^2)),0) // Epanechnikov Kernel*
*gen w=max((1/(2*3.14159265)^(1/2))*exp((-1/2)*((date-mdy(08,20,2005))/60)^2),0) // Gaussian Kernel

gen polya= (date-mdy(02,23,2013))*(1-EL)
gen polyb= (date-mdy(02,23,2013))*EL

*** Hourly 
set more off
capture drop _I*
qui xi i.station*polya i.station*polyb   i.station*i.dow 
local covars EL polya polyb ln_210 ln_gp vis vis2 prec prec2 _I*
forval h=5/9{
	foreach w in w55 w60 w65{
		foreach ln in HV ML{
			local cond type=="`ln'" & `w'!=0 & hour==`h'
			qui{
			 reg ln_flow `covars' [aw=`w'] if `cond', cluster(weekid)
			estimates store flow_mean_b`w'_`ln'
			 reg ln_speed `covars' [aw=`w'] if `cond', cluster(weekid)
			estimates store speed_mean_b`w'_`ln'
			
			
				}
			}
		}
	*cd "\\rschfs1x\userrs\arw227_RS\Documents\expresslanes\Fall2014\tables"
	esttab flow_mean_bw55_HV  flow_mean_bw60_HV flow_mean_bw65_HV flow_mean_bw55_ML  flow_mean_bw60_ML flow_mean_bw65_ML ///
		speed_mean_bw55_HV  speed_mean_bw60_HV speed_mean_bw65_HV speed_mean_bw55_ML  speed_mean_bw60_ML speed_mean_bw65_ML ///
			using .\results\revision1\RDD_results_hrly_mean_`h'.csv, replace cells(b(star fmt(3)) se(par fmt(3)))  ///
			   legend label varlabels(_cons constant)              ///
			   stats(r2 N) keep(EL ln_210 ln_gp vis vis2 prec prec2)
			   
	
	}
	
*** Period
	
	
	set more off
capture drop _I*
qui xi i.station*polya i.station*polyb  i.station*i.dow i.station*i.hour
local covars EL polya polyb ln_210 ln_gp vis vis2 prec prec2 _I*
	foreach w in w55 w60 w65{
		foreach ln in HV ML{
			local cond type=="`ln'" & `w'!=0 & period=="AM"
			qui{
			 reg ln_flow `covars' [aw=`w'] if `cond', cluster(weekid)
			estimates store flow_mean_b`w'_`ln'
			 reg ln_speed `covars' [aw=`w'] if `cond', cluster(weekid)
			estimates store speed_mean_b`w'_`ln'
			
		
				}
			}
		}
 	esttab flow_mean_bw55_HV  flow_mean_bw60_HV flow_mean_bw65_HV flow_mean_bw55_ML  flow_mean_bw60_ML flow_mean_bw65_ML ///
		speed_mean_bw55_HV  speed_mean_bw60_HV speed_mean_bw65_HV speed_mean_bw55_ML  speed_mean_bw60_ML speed_mean_bw65_ML ///
			using .\results\revision1\RDD_results_mean_AM.csv, replace cells(b(star fmt(3)) se(par fmt(3)))  ///
			   legend label varlabels(_cons constant)              ///
			   stats(r2 N) keep(EL ln_210 ln_gp vis vis2 prec prec2)
	
   
   
   
   
